%% Posterior Function

function PS = posterior_fun_MA(xx,ydata,prior,NS,Nx,NC,tr,setmin)

lnL_sbl_xC = 189*ones(NC,1);
xx = xx';

if tr == 1  % --> parameters in R, transform it to right domain
    alphap11 = exp(xx(1))/(1+exp(xx(1))); thch.p11 = alphap11*prior.p11lb + (1-alphap11)*prior.p11ub;
    alphap22 = exp(xx(2))/(1+exp(xx(2))); thch.p22 = alphap22*prior.p22lb + (1-alphap22)*prior.p22ub;
    alphamu1 = exp(xx(3))/(1+exp(xx(3))); thch.mu1 = alphamu1*prior.mu1lb + (1-alphamu1)*prior.mu1ub;
    alphaDmu = exp(xx(4))/(1+exp(xx(4))); thch.Dmu = alphaDmu*prior.Dmulb + (1-alphaDmu)*prior.Dmuub;
    alphasg  = exp(xx(5))/(1+exp(xx(5))); thch.sg  = alphasg *prior.sglb  + (1-alphasg) *prior.sgub;        
    thch.mu2 = thch.mu1 + thch.Dmu;
else
    thch.p11 = xx(1);
    thch.p22 = xx(2);
    thch.mu1 = xx(3);
    thch.Dmu = xx(4);
    thch.sg  = xx(5);        
    thch.mu2 = thch.mu1 + thch.Dmu;
end
thch.Pmat = [thch.p11 ,  1-thch.p11;...
             1-thch.p22,  thch.p22];         
thch.F = [thch.sg, -thch.sg];
thch.R = 0;
thch.Q = 1;
thch.A = [0, 0 ; ...
          1, 0];           
thch.G = [1; 0];
thch.muvec = [thch.mu1, thch.mu2] ;
thch.x00 = zeros(Nx,1);
thch.P00 = eye(Nx);

% run Compute_XI_ss
Compute_XIss
thch.XI_00 = XI_ss;

PP = 1; % PP = prior_fun(thch,prior);

for ii = 1:NC
    yaux = ydata(:,ii);
    in   = find(isnan(yaux) == 0, 1, 'first');
    fn   = find(isnan(yaux) == 0, 1, 'last');
    y    = yaux(in:fn);
    [~, ~, ~,~, ~, lnL_sbl_xC(ii)] = Kim_Filter_MA(thch,y,NS,Nx);    
    lnL_sbl_xC(ii)                 = log(PP) + lnL_sbl_xC(ii) ;    
    
end

PS = sum(lnL_sbl_xC);

if setmin == 1
    PS = - PS;
end

return